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Abstract. The deduction from solar flare X-ray photon spectroscopic data of 
the energy dependent model-independent spectral index is considered as an inverse 
problem. Using the well developed regularization approach we analyze the energy 
dependency of spectral index for a high resolution energy spectrum provided by 
Ramaty High Energy Solar Spectroscopic Imager (RHESSI). The regularization tech- 
in : nique produces much smoother derivatives while avoiding additional errors typical 

of finite differences. It is shown that observations imply a spectral index varying 
significantly with energy, in a way that also varies with time as the flare progresses. 
■ The implications of these findings are discussed in the solar flare context. 
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o 
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1. Introduction 

o 

Hard X-ray spectroscopy is considered to be an important tool for the 
study of high energy processes at the Sun (e.g. Brown, 1971; Lin and 
Hudson, 1976; Aschwanden, 2002), yielding vital, direct information on 
5_i ■ fast electron populations in flares. Specifically, the spatially integrated 

czj , X-ray spectrum may be viewed as the convolution in electron energy 

of the bremsstrahlung cross-section and the mean electron flux in the 
source region (Brown, Emslie and Kontar, 2003). The functional form of 
^ \ the photon energy spectrum then contains information on the energy- 

5^ ' dependence of the mean electron flux. This information may in turn 

be interpreted to bear on pictures of flare electron acceleration and 
propagation. 

Even the earliest observations (e.g. Kane and Anderson, 1970) 
revealed the overall power-law form of the photon spectrum /(e) (pho- 
tons cm -2 s _1 keV^ 1 ) above about 10 - 20 keV: /(e) ~ e~ 7 , for some 
7 > 0. In the simplest possible (isotropic, non-relativistic) treatments 
of bremsstrahlung production, a power-law photon energy spectrum 
implies a source mean electron flux F(E) (electrons s -1 ) which is also 
a power-law in electron energy E. Further, simple assumptions about 
the post-acceleration propagation of electrons in the source (e.g. thin 
target, in which all electrons escape with negligible energy loss, or thick 
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target, in which all electrons stop completely in the source) then lead 
to F(E) ~ E~ s , with simple, linear relationships in each of the thick 
and thin target cases between 5 and 7 (Brown, 1971). 

A pure power-law form of F(E) would be an important clue to the 
nature of the acceleration process (Korchak, 1971). Power-law electron 
energy distributions are found commonly in nature and are a natural 
consequence of certain acceleration mechanisms (Miller et al, 1990; 
Petrosian et al., 1994). The distinctive feature of a power law, however, 
is the absence of any natural scale so it is probably true that any 
deviations from power-law behaviour give more decisive clues to the 
acceleration mechanism than the power law itself. Any such features 
would, of course, be reflected in deviations from power-law behaviour of 
/(e) - see Lin and Schwartz (1987) for an example. In particular, F(E) 
of power-law form only up to some maximum energy E max would result 
in an /(e) which falls off more and steeply as e = E max is approached. 

Even for an accelerated pure power-law F(E), electron trans- 
port and/or radiation physics may produce deviations from power- 
law behaviour in the observed spatially integrated photon spectrum. 
Anisotropy of the mean electron distribution combined with the di- 
rectionality of bremsstrahlung emission for deka-keV electron energies 
would also have this effect, possibly with diagnostic potential for the 
source electron angular distribution (Massone et al, 2004). Compton 
back-scattering of photons from the photosphere ("X-ray albedo") will 
also distort a power-law photon spectrum (Bai and Ramaty, 1978). 

In a thick target, several physical processes may alter the analysis 
that leads to a power-law /(e) from a power-law F(E), e.g. nonuniform 
ionisation (Brown, 1973; Kontar et al, 2003). As emitting electron en- 
ergies approach thermal energies, the growing importance of velocity 
diffusion (as opposed to systematic mean slowing-down) will cause the 
mean electron distribution, and thus the photon spectrum, to deviate 
from a power law (Galloway et al., 2004). Noncollisional energy losses 
could have a wide variety of consequences for the photon spectrum 
(Brown and MacKinnon, 1985). 

It seems clear that deviations from power-law /(e) hold significant 
diagnostic potential. We may define 7(e) = —(e/l)dl/de as an energy- 
dependent spectral index, identically constant for a pure power-law 
/(e), but likely to be informatively non-constant much of the time. 
Most earlier X-ray detectors had insufficient photon energy resolution 
to fully explore this potential, however, but we may now realistically 
aim to calculate 7(e) numerically from RHESSI (Lin et al, 2002) data. 

With the simplest possible approximation (Kramers, 1923) to the 
bremsstrahlung cross-section we see further the relationship between 
the spectral index and the mean electron flux F(E), in particular il- 
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lustrating physically the connection between derivatives and an inverse 
problem. The photon flux is a convolution of the mean electron flux 
and the cross-section (Brown, 1971). Making Kramers' approximation 
we can write 

e 

where Qq is a constant. From the observed flux /(e) we need to obtain 
the mean electron flux. Differentiating both parts of © immediately 
leads us to the following 

F(E)\ E=e = e -^[l(e)-l] (2) 

Since the photon flux is given from observations, Equation © shows 
that the spectral index uniquely determines the mean electron flux, 
at least in the limit that Kramers' approximation applies (Kramers, 
1923). 

Real data always come with noise, however, and numerical differ- 
entiation of data is a noise-amplifying process. Thus there is a need 
to establish reliable procedures for obtaining 7(e) (and more generally, 
dl/de) safely from data. Here we demonstrate that regularisation tech- 
niques (Groetsch, 1984; Hanke and Scherzer, 2001) can be applied for 
our purposes and yield such a robust procedure. The resulting smooth 
function may then be safely differentiated to obtain a best estimate of 
7(e). In Section 2 we give the formal demonstration of such an approach, 
in the process clarifying the sense in which this gives a "best" estimate. 
Section 3 applies this technique to obtain 7(e) from RHESSI data for 
the flare of 26 February, 2002 and study its time-dependence. Section 4 
gives conclusions and discusses physical implications of the work carried 
out here. 



2. Derivative as an inverse problem 

Let us assume that we have a smooth function y(x) over the interval 
%oi < x < xo2- We have a finite sample yi of measured values of this 
function, obtained over some grid 201 = xq < x\ < ... < Xi < ... < 
x n = xq2 with mesh size Ax. The noisy data set has an error 

\Vi -y{xi)\ < 5y (3) 
where 5y is an uncertainty of measurement. 
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We want to find the best smooth estimate of the derivative y'(x) 
using the given data set V x E [xoi,xo2]- The two point finite differ- 
ence estimate is readily available from the Cauchy expansion with the 
following bound 

< 0{Ax + 5y/Ax), (4) 

where the first and second terms in the right hand side represent con- 
sistency and propagation errors respectively (Groetsch, 1984). The first 
term in the right hand side comes from the discreteness of the data set, 
while the second is connected with the errors of the data 5y. 

Equation @ explicitly shows that for decreasing Ax the error in 
the estimated derivative deteriorates rapidly, whereas for increasing 
Ax the error grows only linearly. The propagation error can be sub- 
stantially suppressed by taking larger Ax, a procedure also known as 
regularization by coarse discretization. The extreme case of discretiza- 
tion corresponds to a linear fit. Obviously, this is far from desirable: 
although we have minimised the effects of error propagation we have 
only a single estimate of y'(x) across the whole of our range of Xj. On 
the contrary, we want to extract as much information as possible. The 
finite difference approach also leads to a non-smooth estimate of y'(x) 
- i.e. the resulting derivative is piecewise continuous. Note that the 
right-hand side of the Equation @ reaches its minimum value which 
is 0{-\f5y) when Ax ~ \ftyy~. 

Let us now assume the existence of a function /(x), which is a close 
approximation to the data set: 

-^£(y 4 -/(^) 2 <(5y) 2 , (5) 

1=1 

Below we describe a detailed procedure for the construction of /(x). 
We have to suplement this procedure with boundary conditions, the 
values of f(x) at two values of x. The values yi at any two of the 
values Xi would suffice in principle. It might under some particular 
circumstances (e.g. suspected outliers) be appropriate to choose two of 
the interior points, but in what follows we simply choose the values 
at the endpoints, so that boundary conditions become: /(xoi) = yo, 
and /(X02) = y n - Let us also assume that /(x) is a square integrable 
function over interval x G [xoi,xq2], so we can introduce a norm 




Hanke and Scherzer (2001) describe an optimal approach to con- 
structing a smooth estimate of the derivative y'(x). Following them, 



Vi+l ~ Vi 



Ax 



y'(xi) 
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we also require that the function f(x) has a smooth derivative, in the 
sense that 

||/"||=min (7) 

This means that we seek the function f(x) whose second derivative 
has the smallest norm. Roughly speaking, if we represent f(x) locally 
by the first two terms of its Taylor expansion, f"(x) gives an estimate 
of the error on f'(x) and minimising ||/"|| minimises the error in the 
resulting estimate for /'. 

More specifically, let us estimate the error of our approximation if 
our basic conditions ©,0 are met. The propagation error for the 
derivative /' can be estimated by considering the following norm 

\\f'-i/\\ 2 <\\f-v\\ II/" -/II (8) 

where we integrated by parts and used the Cauchy-Schwarz inequality. 
The first factor in the right hand side is bounded by our requirement 
The second factor is also limited as soon as ||/"|| < ||y"||> which 
is true as long as equation © is satisfied. Now using the Minkowski 
inequality we arrive at the following result 

\f'(x i )-y'(x i )\<0(Ax+^) (9) 

where we added consistency error in the same manner as in Q. 

Comparing @ and (JU we see the drastic difference. The former 
shows no growth for small values of Ax: fitting a suitably smooth 
function prior to estimating the derivative eliminates the uncertainty 
associated with estimates based on discrete bins. Moreover, the error 
© is bounded by the minimum in (flf). 

The result of this formal exercise can be formulated in terms of opti- 
mization (Hanke and Scherzer, 1998). Collecting our two requirements 
©,0 into a single equation, we require to minimise the functional: 

$(/)^^X>- r f'(0^-yo) 2 + X\\f"(x)\\ 2 (10) 

n 1 i=1 Jx 

among all smooth functions / with /(eoi) = y(^oi)> f( x 02) = 2/(^02), 
where A is so that 

-^-f> - h{xi)f = (5y) 2 (11) 

n — 1 r— i 

i=i 

where f\(xi) is a solution of minimum problem (|10j) . This is formally 
equivalent to an inverse problem and the solution method employed 
is well-known as Tikhonov regularization (Tikhonov, 1963). The equa- 
tion (|llj) is also known as a discrepancy principle. It is interesting 
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to compare the problem (|10[) with the inversion of the photon data 
using second order regularization (see Kontar et al, 2004). The two 
problems are virtually identical. The only difference is the operator: the 
bremsstrahlung cross-section used in case of inversion of photon data, 
while our operator is an integrator. For the numerical implementation of 
this exercise, it is trivial to write the formal solution of the problem 
using Generalized Singular Value Decomposition of our two operators: 
integrator and second order derivative (Kontar et al, 2004). 

The solution of the problem (|10lll|) presents a regularized solution of 
the problem for f(x). Figure^shows an example of this process applied 
to the function y(x) = cos(x), with a perturbation added in the form of 
random noise at the level of 1% of the value of y. The upper panel shows 
the actual data points generated in this way, together with the original 
function to emphasise how close the "data" points generated actually 
are. The lower panel shows values of the derivative calculated just by 
simple differencing of the data points, and the smooth curve generated 
by first regularising the data in the way described above. Judging from 
this, the derivative resulting from the regularization procedure is clearly 
much better. Compared to this, simple differencing of data with even 
this modest level of noise destructively affects the result. 



3. Energy variation of spectral index in a solar flare 

We defined the energy-dependent spectral index as a logarithmic 
derivative: 7(e) = — dm //dine, where /(e) is an observed spectrum 
given as a data set Ij for every e,. Choosing quantities cc, and y% as 
follows 

ln(Ji) -> yu ln(ei) -» x, (12) 

we find our desired spectral index to be 7(e) = —f'(x), where fix) is 
determined from the measured values yi as described in the previous 
section. We note here that since the bremsstrahlung photon spectrum 
is a convolution of cross-section and electron flux (Brown, Emslie and 
Kontar, 2003), logarithmic derivative of photon spectrum should always 
exist. 

We now give a first, illustrative application using the RHESSI data 
from the GOES M-class flare that happened on February 26, 2002 
around 10:27 UT (Figure EJ). This event provides good photon count 
statistics while below the level of nonlinear features of the instrumental 
response such as pulse pile-up in the detectors (Smith et al, 2002). 
(These features are not fully understood and their details lie outside 
the scope of this paper.) Using standard software tools in SPEX we 
extracted the photon spectral flux (photons keV -1 s _1 cm~ 2 ) in the 
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Figure 1. Numerical example. Function I(x) = cos(a;) (solid line) and data points 
(crosses) (upper panel) and its derivative (lower panel). Solid line is a regularized 
derivative for data points. The dash line shows the actual values of finite differences. 

range above lOkeV accumulated in seven front segments out of nine 
detectors. We omitted detectors 2 and 7 due to low operational energy 
resolution at the moment of the observation. 

Figure 01 shows the spectrum of the solar flare observed by RHESSI 
on February 26, 2002 for the time interval 10:26:40-10:27:00 UT near 
at the peak of the flare. The lower panel shows the spectral index, 
7(e), calculated via the regularisation process described above. To con- 
struct error bars we first generated two further sets of "data", one 
la above and one la below the original data points I fa), where a is 
the uncertainty on the photon spectrum. Then we applied the same 
regularisation procedure to these two new data sets. The resulting con- 
fidence interval in 7(e) allows us to highlight several features of interest, 
particularly in comparison with a other parametric descriptions of the 
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Figure 2. Temporal variation (4 seconds cadence) of the count rates in seven front 
RHESSI segments for the February 26, 2002 solar flare. The vertical lines show five 
20 second selected accumulation intervals for spectral analysis. 



spectrum and by implication its derivatives, for example an isothermal 
plus broken power-law spectral fit (Fig. |SJ) - cf. Holman et al. (2003). 

Most immediately obvious (Fig. is the nonconstancy of the spec- 
tral index within the energy range studied (10-100keV). At the lowest 
photon energies the regularized spectral index has a higher value than 
the index suggested by the isothermal plus power-law fit. This be- 
haviour of spectral index may indicate X-ray emission from plasma 
with a range of temperatures. Similar suggestions of non-isothermality 
have been made by Piana et al. (2003), analyzing the mean electron 
flux spectrum deduced in the July 23 2002 flare and comparing to an 
isothermal and broken power-law fit (Holman et al, 2003). The regu- 
larized spectral index 7(e) also departs significantly from the broken 
power-law fit near the break energy 49 keV and above 70 keV. At 
the highest energies studied, 7(e) approaches 4, larger than the value 
found assuming a broken power-law. Above the thermal component, 
the spectral index grows approximately linearly with energy. 

Figure [5] shows that for all time intervals the regularized spectral 
index shows a minimum, the position of which grows with time. There 
is a clear minimum for the first time interval at 10:26:00 UT near 17 
keV, while at 10:27:20 UT the minimum is less clear and as high as 
50 keV. Additionally, the spectral index varies less with energy as the 
flare progresses. The spectral index growths with energy by as much 
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as 3.0 for the first time interval, but it varies only by 0.1 for the last 
interval. 

Earlier data from scintillation detectors established a relationship 
between spectral hardness (parametrised via a single value of 7 derived 
from fitting across the available energy range) and total photon flux 
above some energy (Kane and Anderson, 1970). Specifically, 7 obtained 
in this way is anti-correlated with total photon flux, so that the flare 
X-rays, viewed crudely, have the hardest spectrum when they are most 
intense. This "soft-hard-soft" spectral behaviour has been confirmed 
with RHESSI using ratios of photon flux in fairly broad bands (Fletcher 
and Hudson, 2002; Hudson and Farnik, 2002). Figure El shows the full, 
complex behaviour behind this simplified description. Values of 7(e) are 
indeed lowest, for any given e, around the peak of the event, but the 
decrease before and increase after this time do not occur at the same 
rate for all e. The flare starts with spectral index strongly dependent 
on energy, while a single power-law fit produces a mean spectral index 
which is a combination of small spectral index around 19 keV and 
relatively large 7(e) at 70 keV. In fact the minimum of the spectral 
index shows a tendency to grow as the flare progresses. Part of the 
transition from soft to hard spectra can thus be understood in terms of 
a reduction in the degree of variation with energy of 7(e) in the range 
20 - 100 keV. 



4. Discussion 

The illustration above shows vividly how far real spectra depart from 
the idealised power-law form. In this section we discuss some possible 
physical implications. 

In standard models of hard X-ray production (e.g. Aschwanden, 
2002) the X-ray emitting region is located well above the dense pho- 
tosphere. X-rays emitted downwards will scatter on electrons of the 
dense layers of the solar atmosphere, finding their way back to an 
observer. The scattered flux, a function of the flare location and primary 
spectrum, may produce a detectable alteration of the net spectrum in 
the range 15-60 keV (Bai and Ramaty, 1978). Assuming the primary 
photon spectrum to be a simple power law (Bai and Ramaty, 1978) the 
photospheric albedo will peak at energies in the 20 - 40 keV range, with 
the reflected flux ~ 30 — 70% of the primary flux. As a result, even if 
the primary spectrum is a pure power-law, the observed spectrum is 
flatter before the albedo peak energy and steeper above. This should be 
seen in the energy variation of the spectral index. For the February 26, 
2002 flare, at heliocentric angle ~ 75°) and with photon spectral index 
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Figure 3. Spectrum of February 26, 2002 solar flare for 10:26:40-10:27:00 UT (upper 
panel crosses). The low panel shows the inferred spectral index with errors. The thick 
solid line presents the spectral index based on the isothermal plus broken power-law 
fit with 7; om = 2.8,7u P per = 3.5 with the break energy 49 keV. 



s=! 3, the results of Bai and Ramaty (1978) suggest that the albedo 
contribution to the observed spectrum will maximise near 35 keV, 
implying a lower spectral index below this energy and higher above. 
This tendency can indeed be seen in Figure ©• Nevertheless, albedo 
alone cannot account for the temporal variations shown in Figure 
however - these must reflect properties of the primary X-ray spectrum, 
and thus of the emitting electron population. 

Suppose that there is a high energy cut-off (or substantial soften- 
ing) in the bremsstrahlung-emitting electron spectrum, at some energy 
Emax- For energies close to E max the spectral index will be a growing 
function of e. Since the photon flux at a given energy is an integral 
over all electrons above this energy, the change of the spectral index 
due to a high energy cut-off can be seen at photon energies e well below 
E m ax (Kontar et al, 2004). During the first time interval the variation 
with e of the spectral index is most significant, while the photon flux 
above 50 keV is very low (Figure [2J). These facts together suggest a 
small number of electrons above 50 keV at this early time, while later 
in the flare the number of high energy electrons (above 100 keV) grows 
substantially. This is at the very least a plausible interpretation of 
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Figure 4- The residuals and the accumulated residuals for the time interval 
10:26:40-10:27:00 UT (Fig. 0J. The dash lines in the low panel present 3a level 
assuming residuals are statistically independent. 
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Figure 5. Temporal variation of energy dependent spectral index. Each line corre- 
sponds to one time interval. 
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the findings of Figure EJ suggesting that the "soft-hard-soft" pattern of 
spectral evolution really reflects changes in the highest electron energies 
present at different times. 

A number of other physical processes may produce detectable varia- 
tion in spectral index, even given an injected power-law electron energy 
distribution. For example, it was suggested (Brown, 1973; Kontar, 
Brown, and McArthur, 2002) that observed deviations from power-law 
can be interpreted as a manifestation of nonuniform target ionization. 
However, the possible variations of spectral index are unlikely to exceed 
0.6 (Brown, Emslie, and Kontar, 2003). Moreover, the minimum of 
the spectral index if associated with chromospheric depth should show 
modest growth after the impulsive phase of the flare (Kontar et al, 
2003), while the minimum of the spectral index shows no correlation 
with X-ray flux. 



5. Conclusions 

We have shown how best to estimate numerically the energy-dependent 
X-ray spectral index 7(e) from a set of observed photon fluxes at dis- 
crete energies, in a way that avoids the noise amplification of numerical 
differentiation. Formally the necessary procedure is equivalent to regu- 
larization. In consequence we have shown that the spectral index of at 
least one solar flare can have large variations with energy. The spectral 
index shows a clear minimum in the range 17-50 keV and the value 
of this minimum decreases as the flare progresses. The spectral index 
also tends towards energy-independence as the flare progresses, sug- 
gesting that something more like an ideal power-law photon spectrum 
is eventually attained. The origin of such variation can be connected 
with variability in the highest energy of the X-ray producing electrons. 
For later times the high energy cut-off is higher, thus producing a more 
uniform spectral index. This seems a simple, plausible explanation of 
the spectral behaviour found here; in particular it serves to illustrate 
the value of obtaining and studying 7(e). 
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